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We study the magnetization for the classical antiferromagnetic Ising model on the Shastry- 
Sutherland lattice using the tensor renormalization group approach. With this method, one can 
probe large spin systems with little finite-size effect. For a range of temperature and coupling con- 
stant, a single magnetization plateau at one third of the saturation value is found. We investigate 
the dependence of the plateau width on temperature and on the strength of magnetic frustration. 
Furthermore, the spin configuration of the plateau state at zero temperature is determined. 

PACS numbers: 75.60.Ej, 05.50.+q, 75.10.Hk, 05.10.Cc 



I. INTRODUCTION 



The frustrated spin systems have attracted much at- 
tention over last decades since very rich physics can ap- 
pear in these systems.— Some interest in such systems 
is concentrated on fascinating sequence of magnetiza- 
tion plateaus at fractional values of the saturation mag- 
netization, which was first observed in two-dimensional 
spin-gap material SrCu2 (603)2'^ This compound can be 
described well by spin-1/2 antiferromagnetic Heisenberg 
model on the frustrated Shastry-Sutherland lattice (or 
the orthogonal-dimer lattice))^ as shown in Fig. [TJ Be- 
sides the previously discovered plateaus at 1/3, 1/4 and 
1/8 of the saturated magnetization, evidence in favor of 
more fractional magnetization plateaus down to values as 
small as 1/9 has been reported recentlyj^i^ Stimulated 
by the discovery of magnetization plateaus, various theo- 
retical and experimental explorations have been devoted 
to the properties of the Shastry-Sutherland model in a 
magnetic field.— i^i^iiiS 

Similar phenomena of magnetization plateaus is also 
observed in rare-earth tetraborides RB4. The mag- 
netic ions of these compounds are again located on a 
lattice that is topologically equivalent to the Shastry- 
Sutherland lattice . ^^i^^'^'^i^'^'^^i^^'^^ In particular, magne- 
tization plateaus at small fractional values (1/7, 1/9 .. . 
of the saturation magnetization) are reported in the com- 
pound TmB4.-'^^-^^ Because fully polarized state can be 
reached for experimentally accessible magnetic fields, this 
compound allows exploration of its complete magneti- 
zation process. Note that, due to large total magnetic 
moments of the magnetic ions, this compound can be 
considered as a classical system. Moreover, because of 
strong crystal field effects, the effective spin model for 
TmB4 has been suggested to be described by the spin-1/2 
Shastry-Sutherland model under strong Ising (or easy- 
axis) anisotropy.^^ Thus, studying the Ising limit is the 
first step toward a complete understanding of the mag- 
netization process for this material. 

In the presence of a finite magnetic field h, the to- 
tal energy of the antiferromagnetic Ising model on the 
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FIG. 1: The Shastry-Sutherland lattice. J bonds (dashed 
lines) are the exchange couplings along the edges of the 
squares and J' bonds (solid lines) are the diagonal intra-dimer 
couplings. 



Shastry-Sutherland lattice is given by 



(1) 



with exchange couphngs J, J' > 0. Here, Si — ±1/2 
denotes the z-component of a spin-1/2 degree of freedom 
on site i of the square lattice. The first sum extends over 
all nearest neighbor bonds, and the second sum runs over 
next-nearest neighbor bonds in every second square, as 
indicated in Fig. [1] Even for this simplified case, differ- 
ent conclusions for the magnetization curve have been 
reached. In Ref. [13, a single plateau at 1/2 of the sat- 
uration magnetization is found based on analyzing a fi- 
nite system with 16 spins only. However, when larger 
system sizes up to 18 x 18 spins are considered, a dis- 
tinct plateau at 1/3 of the saturation magnetization is 
obtainedii^ The discrepancy may come from the effect of 
finite lattice sizes. As noted by the authors of Ref. [H, for 
finite systems, inappropriate lattice sizes and boundary 
conditions can frustrate certain magnetization patterns, 
and hence lead to rather different magnetization curves 
which do not correctly represent the behavior in the ther- 
modynamic limit. For example, the plateau at 1/3 of the 
saturation magnetization is not allowed for systems of 



4x4 and 8x8 spins, even though it does describe the 
true magnetization process for systems in the thermody- 
namic hmit. 

In order to check theoretically if other reported mag- 
netization plateaus at small fractional values can be sta- 
bilized in the current model, unbiased large-scale calcula- 
tions are called for. This is because the unit cells of mag- 
netization profiles inside high-commensurability plateaus 
are usually quite large, calculations for systems of fi- 
nite sizes may prevent reliable predictions for these cases. 
Therefore, to avoid the frustration for certain magnetiza- 
tion plateaus coming from geometric constraints, and in 
particular to uncover the possibility of plateaus at small 
fractional values, analyzing systems of large enough sizes 
are necessary. 

Lately, based on ideas from quantum information the- 
ory, the tensor renormalization group (TRG) method is 
developed," which can efficiently calculate quantities of 
classical systems of very large sizes. This technique can 
in principle be applied to any classical lattice with local 
interactions as long as the partition function can be ex- 
pressed as a tensor network.^" Because the accuracy can 
be systematically improved by increasing the cutoff on 
the index range of the tensors, highly precise quantities 
can be calculated under the TRG approach even in the 
thermodynamic limitii^i^ii^ Therefore, the TRG method 
is one of the most suitable ways to study the magneti- 
zation process of the classical frustrated spin systems in 
the thermodynamical limit. 

In the present work, the magnetization process of the 
spin-1/2 Shastry-Sutherland model in the Ising limit is 
investigated by employing the TRG approach J^s^i^^ We 
find that the magnetization curve exhibits exactly one 
plateau at 1/3 of the saturation value. Our results are 
in accordance with the findings in Ref. Furthermore, 
phase diagrams in the (/i, T) plane for a typical magnetic 
coupling ratio J' /J = 1 and in the {h, J') plane for a 
particular temperature T/J = 0.2 are obtained. Since 
there is no evidence for the presence of any additional 
plateaus for the spin-1/2 Shastry-Sutherland model in 
the Ising limit, to explain the experimental results, one 
must go beyond this simple model. 

This paper is organized as follows. In Sec. II, the TRG 
approach is outlined briefly. In Sec. Ill, we apply this 
method to investigate the magnetization process of the 
Shastry-Sutherland model in the Ising limit. The spin 
configuration of the plateau state at zero temperature is 
discussed in Sec. IV. Sec. V is the conclusion. 



II. TRG FORMULATION 

Before applying the TRG method of Levin and NavCfi^ 
we first explain how to express the partition function 
of the present model as a tensor network. One possi- 
ble way is to rewrite the total energy in Eq. ([T]) as a 
summation over the energies of plaquettes with diagonal 
bonds. The energy of, say, the plaquette A with spins 
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FIG. 2: (Color online) Checkerboard decomposition of the 
Shastry-Sutherland lattice and the corresponding tensor net- 
work. 



si, S2, S3, S4 on its corners is given by (see Fig. [2]) 

e^(si,S2, 33,54) = J{S1S2+ S2S3 + S3S4 + S4S1) 

+ J'S2S4 - -^(Sl + S2+ S3+ S4) .(2) 

The rank-four tensors are defined as the Boltzmann 
weights for these plaquettes. For example, 

T'ai,a2,a3,a4 = exp[-/3e^ (si , S2 , S3 , S4)] (3) 

with f3 being the inverse temperature and the indices 
ai = Si + 3/2 running over 1 and 2. Afterwards, the 
partition function can be rewritten as a sum of tensor 
products in the following way, 

Z = ^e-'^^^^^'" 

{S^} 

= tTr(T^T^---) . (4) 

Here the tensor trace (tTr) means that all indices on the 
connected links in the tensor products are summed over. 
As a result, the partition function of the Ising model on 
the Shastry-Sutherland lattice is transformed to a tensor 
network as shown on the right-hand side of Fig. [2] 

As discussed in Refs. [l9ll2lll2^ the tensor network can 
be coarse-grained in an iterative fashion to reduce the 
load of computation. At the mean time, the accuracy 
can be controlled by a parameter of cutoff D^ut ■ Here we 
outline the process briefly. Each step of the renormaliza- 
tion consists of two operations: rewiring and decimation. 
After one step of the renormalization, the number of sites 
in the tensor network is reduced by half. Eventually, the 
system is reduced to only four sites (four T's) and the 
partition function can be calculated with ease. 

Rewiring - By viewing the rank-four tensor as a ma- 
trix, say Mi^a2,a3)Xa4,ai) = ^^^i, as, 03, ^ud with the 
help of singular value decomposition (SVD), Af — UAV\ 
the rank-four tensor can be decomposed to two rank- 
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FIG. 3: (a) Rewiring: the original rank-four tensors are de- 
composed to two rank-three tensors, (b) Decimation: the new 
tensor T' is obtained by summing over the indices around the 
square. 



three tensors. That is [see Fig.[2Ia)], 
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(similarly for S''^ and S^), in which A.^ are 
the singular values, and U, V are the unitary matrices in 
SVD. If each index of the original rank- four tensor has 
D possible values, then there should be terms in the 
summation of Eq. ([5]). In practice, the tensor is approx- 
imated by keeping only the largest Dcut singular values 
and the corresponding singular vectors. Apparently, the 
cutoff needs to be chosen such that the result converges 
with little Dcut-dependence. 

Decimation - After rewiring, the dashed lines in 
Fig. [Slja) can be closed to build a new rank-four ten- 
sor, T' [see Fig. [3Ib)]. This is achieved by the following 
operation, 

(6) 



T' 

71,72,73,74 



V-^72-'7l 74-'73/ 



where the square matrices (SM = Sf^ . After 

V / 0,0' ^ ' 

such a contraction, one obtains a new tensor network that 
is half of the size (see Fig. HJ. Afterwards, the renormal- 
ization can be carried out iteratively until there are only 
four sites left. 
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FIG. 4: Under the TRG procedure, a tensor network is trans- 
formed into a coarse-grained tensor network. 



We note that, to prevent the computation from diverg- 
ing, one needs to normalize the new rank-four tensor at 
each step of the renormaliztion. At the beginning, we 
factor out the maximal value Wq of the tensor elements 
of j'A/B ^ Tq^^ to obtain a normalized tensor T^^^ . 
After the first step of the renormalization-group (RG) 
transformation on T^^^, a renormalized tensor T' = Ti 
is reached. Now we choose the normalization factor to 
be 1^1 = X:^ax^^ax such that Ti = Wifi, where X:^^^ 
and X^cix ^I'S the largest singular values of the two de- 
compositions in Eq. ([5]). 

The factorization and RG transformation arc then 
iterated, so that at the nth step we have a tensor 
Tfi = WnTn- Thus, for the Shastry-Sutherland lattice 
of = 2"+'^ sites (and with N/2 tensors in the original 
tensor network) , after n steps of the RG transformation, 
one has 
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Since the last tensor-trace term in Eq. ([7]) remains finite, 
its contribution to the free energy can be neglected for 
a large enough system. The free energy per site thus 
becomes 
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(8) 



Once the free energy is obtained, one can proceed to cal- 
culate the magnetization. The results are shown and 
discussed in the following sections. 



III. NUMERICAL RESULTS 

In this section, we present the numerical results on 
the magnetization plateau and related phase diagrams. 
Throughout the region being explored, we find only one 
magnetization plateau at m/ms = 1/3, where m denotes 
the magnetization and its saturation value. Unless 
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otherwise mentioned, the size of the system is 2 
with periodic boundary condition. That is, the number 
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FIG. 5: Magnetization curves for three different tempera- 
tures. The size of the system is 2^" x 2^^". The parameters 
are J' = 1 and Dcut ~ 18. 

of steps of the RG transformation in Eqs. ^ and ([8]) is 
n = 17. The temperature T and the strength of magnetic 
frustration J' are measured in units of J. 

Fig. [5] is a typical diagram of the magnetization curves 
for J' = 1. The curves for three temperatures (T = 0.05, 
0.1 and 0.15) are shown. The size of the system is 
2io ^ and the cutoff Dcut — 18. Current result con- 
verges well against further increase of the system size and 
the cutoff. For example, for T = 0.05, a larger system 
with 2^^ X 2^^ {Dcut = 18) yields a result agrees to the 
sixth decimal place for the most part of the curve. A 
larger cutoff D^ut = 24 (system size 2^° x 2^") shows 
similar accuracy. The result is slightly less accurate near 
the edges of the magnetization plateau but still shows 
no visible difference from the T = 0.05 curve in Fig. [5] 
Compared to other methods, the TRG method is both 
accurate and efhcient for very large systems. 

A more complete scan of the temperature can be found 
in Fig. \6la). Over the whole range of calculation, there 
is only one plateau at m/ms = 1/3. Its width gradually 
shrinks to zero near temperature T = 0.18. The cor- 
responding phase diagram for the 1/3-plateau is shown 
in Fig. [6l^b). The extent of the plateau is determined by 
the locations of maximum slope near its edges, which will 
be denoted as hc,i and hc,2 for the lower and the higher 
critical fields respectively. In Fig. [Hb), we have added 
the theoretical critical fields (/ic,i, ^c,2) = (1, 5/2) at zero 
temperature (details later) . One can see that the numer- 
ical result does extrapolate to the theoretical values as 
temperature decreases. 

In Fig.[7Ka), we show another scan of the magnetization 
with respect to J' and ft, at T = 0.2. At this temperature, 
there is no plateau for small frustration. The plateau ap- 
pears when J' is slightly larger than 1. One can see that 
the widths of the plateaus remain roughly the same for 
J' > 2. Their positions appear to shift linearly with 
respect to the strength of the frustration J'. One can 
see this clearly in the phase diagram of Fig. [7Kb). The 
plateaus are again determined by the locations of max- 



(a) 




0.5 - 

I 1 1 1 

0.05 0.1 0.15 0.2 

T 

FIG. 6: (a) Magnetization versus temperature T and mag- 
netic field h. The parameters are J' = 1 and Dcut = 18. (b) 
Phase diagram of the magnetization plateau. The theoretical 
values of the critical fields at zero temperature are denoted 
by filled circles. 



imum slope. The characters of this phase diagram at 
finite temperature are inherited from its counterpart at 
zero temperature (details later), which is also plotted in 
Fig. [TJb) for comparison. The plateaus at zero tempera- 
ture indeed exhibit a constant width at large frustration 
and a linear shift of the plateau position. Such a behavior 
will be explained in the next section. 



IV. MAGNETIZATION PLATEAU AT ZERO 
TEMPERATURE 

When the temperature is zero, the system is in the 
ground state. If the spin configuration of the ground 
state is known, then the Ising energy of the system can 
be calculated analytically. Afterwards, by comparing the 
ground state energies at different parameters, one can de- 
termine the phase boundaries in the parameter space. In 
this section, we will consider three regimes of magnetiza- 
tion: the unmagnetized state (to = 0), the state of the 
1/3-plateau, and the fully-magnetized state (to/tos = 1). 
It will be shown that the phase boundaries being deter- 
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FIG. 7: (a) Magnetization versus frustration J' and magnetic 
field h. The parameters are T = 0.2 and Dcut ~ 18. (b) 
Phase diagram of the magnetization plateau. Dashed lines 
are the theoretical phase boundaries at zero temperature. 



mined are consistent with the numerical resuhs reported 
in Sec. |TTT1 

In the unmagnetized state with m = 0, we assume that 
the system is either in the Neel state or in the collinear 
state, depending on the strength of the frustration J'. 
These states should be stable when the applied field h is 
small. When the system is in the Neel state [Fig. [H^a)], 
for a unit cell formed by four plaquettes (bounded by 
dashed-dotted lines), there are two sites with spin up and 
two sites with spin down. The nearest-neighbor spins are 
all anti-parallel but the spins connected by the J'-bond 
are parallel. It is not difficult to see that the energy per 
site, including the Zeeman energy (zero here), is, 

1 J' 

in which J = 1. 

For large frustration, the system is more likely to be 
in the collinear state [Fig. [8l^b)]. Again there are two up 
spins and two down spins in a unit cell of four plaquettes. 
Now the energy per site becomes 

im=0 = . (10) 



FIG. 8: (Color online) (a) Spin configuration for the Neel 
state, (b) Spin configuration for the collinear state, (c) Spin 
configuration for the magnetization plateau at m/nis = 1/3. 
The squares with dashed-dotted lines indicate possible choices 
of unit cells. 

By comparing the energies in Eqs. ^ and (fTU)) . one can 
see that the energy of the Neel state is lower (higher) 
than the collinear state when J' < 2 (J' > 2). 

When the applied field increases, the system can un- 
dergo a phase transition to a 1/3-plateau state. There are 
several possible candidates for such a state. In Fig. [5]Jc), 
we show the spin configuration of a state with the lowest 
possible energy. With careful analysis, one obtains the 
following spin energy per site. 



When the applied field is sufficiently strong, then, ir- 
respective of the value of J', the system should be fully 
magnetized. In such a case, it is relatively easy to deter- 
mine the spin energy per site, 

I J' h 

e^-i = -H . (12 

-=" 2 8 2 ^ ' 

By comparing em=o and em/ms=i/3j one can determine 
the boundary between the Neel state and the plateau 
state when J' < 2. The lower critical field hc,i is found 
to be 

h^ i = 2 - J' . (13) 

Similarly, by comparing em=o and em/ms=i/3^ o^i^ has 
the boundary between the collinear state and the plateau 
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state when J' > 2, 

^c,i = -l + y- (14) 

These two straight hnes are indicated as the lower phase 
boundaries at zero temperature in Fig. [Tl^b). 

On the other hand, the upper critical field /ic,2 is ob- 
tained by comparing the energies of the plateau state 
(em/m,=i/3) and the fully magnetized state (em/m^=i), 

/lc,2 = 2+y. (15) 

Such a straight line is also shown in Fig. [Hb) . The area 
bounded by these critical magnetic fields should be the 
maximum width of the plateau when the temperature of 
the system drops to zero. For example, when J' ~ 1, the 
plateau is bounded by (/ic,i, /ic,2) = (1,5/2) at T = 0. 
This agrees nicely with the extrapolation in Fig. [6jb). 

V. CONCLUSION 

The TRG method is applied to explore the plateau in 
the magnetization process for the classical Ising model 



on the Shastry-Sutherland lattice. Systems as large as 
2io ^ sites can be routinely studied with relative ease. 
Therefore, the complications from the finite-size effect 
and its related geometric frustration can essentially be 
avoided. We found a single plateau at m/ms = 1/3 that 
is robust over certain ranges of temperature and mag- 
netic frustration, consistent with the result in Ref. Ts' for 
smaller systems and higher temperatures. The model un- 
der investigation is relevant to the compound TmB4fii 
which is found to have a sequence of plateaus down to 
small fractional values.— We note that the antiferro- 
magnetic transverse exchanges have not been taken into 
account in the current classical model. Therefore, the 
quantum effect caused by these couplings may be essen- 
tial in a full explanation of the observed plateaus. 



Acknowledgments 

M.C.C. thanks the support from the National Science 
Council of Taiwan under Contract No. NSC 96-21 12-M- 
003-010-MY3. MFY acknowledges the support by the 
National Science Council of Taiwan under NSC 96-2112- 
M-029-004-MY3. 



* Electronic address: |mfyang@thu.edu.tw| 

^ Frustrated Spin Systems, edited by H. T. Diep (World Sci- 
entific, Singapore, 2004). 

^ H. Kageyama, K. Yoshimura, R. Stern, N. Mushnikov, K. 
Onizuka, M. Kato, K. Kosuge, C. Slichter, T. Goto, and 
Y. Ueda, Phys. Rev. Lett. 82, 3168 (1999); K. Onisuka, H. 
Kageyama, Y. Narumi, K. Kindo, Y. Ueda, and T. Goto, 
J. Phys. Soc. Jpn. 69, 1016 (2000); H. Kageyama, Y. Ueda, 
Y. Narumi, K. Kindo, M. Kosaka, and Y. Uwatoko, Prog. 
Theor. Phys. Suppl. 145, 17 (2002); K. Kodama, M. Taki- 
gawa, M. Horvatic, C. Berthier, H. Kageyama, Y. Ueda, S. 
Miyahara, F. Becca and F. Mila, Science 298, 395 (2002). 

3 B. S. Shastry and B. Sutherland, Physica B 108, 1069 
(1981). 

* S. E. Sebastian, N. Harrison, P. Sengupta, C. D. Batista, 
S. Francoual, E. Palm, T. Murphy, N. Marcano, H. A. 
Dabkowska, and B. D. Gaulin, arXiv:0707.2075 

^ F. Levy, L Sheikin, C. Berthier, M. Horvatic, M. Takigawa, 
H. Kageyama, T. Waki, and Y. Ueda, Europhys. Lett. 81, 
67004 (2008). 

" M. Takigawa, S. Matsubara, M. Horvatic, C. Berthier, H. 
Kageyama, and Y. Ueda, Phys. Rev. Lett. 101, 037202 
(2008). 

For a review of earlier works, see S. Miyahara and K. Ueda, 
J. Phys.: Condens. Matter 15, R327 (2003), and references 
therein. 

® J. Dorier, K. P. Schmidt, and F. Mila. larXiv:0806.3406:! to 
appear in Phys. Rev. Lett. (2008). 

K. P. Schmidt, J. Dorier, and F. Mila, arXiv:08l0.l596 
° A. Abendschein and S. Capponi, arXiv:0807.1071 , to ap- 
pear in Phys. Rev. Lett. (2008) 
^ S. Yoshii, T. Yamamoto, M. Hagiwara, T. Takeuchi, A. 



Shigekawa, S. Michimura, F. Iga, T. Takabatake, and K. 
Kindo, J. Magn. Magn. Mat. 310, 1282 (2007). 
S. Yoshii, T. Yamamoto, M. Hagiwara, S. Michimura, A. 
Shigekawa, F. Iga, T. Takabatake, and K. Kindo, Phys. 
Rev. Lett. 101, 087202 (2008). 

S. Michimura, A. Shigekawa, F. Iga, M. Sera, T. Taka- 
batake, K. Ohoyama, and Y. Okabe, Physica B 378-380, 
596 (2006). 

S. Yoshii, T. Yamamoto, M. Hagiwara, A. Shigekawa, S. 
Michimura, F. Iga, T. Takabatake, and K. Kindo, J. Phys.: 
Conf. Series 51, 59 (2006). 

F. Iga, A. Shigekawa, Y. Hasegawa, S. Michimura, T. Tak- 
abatake, S. Yoshii, T. Yamamoto, M. Hagiwara, and K. 
Kindo, J. Magn. Magn. Mat. 310, e443 (2007). 
S. Gabani, S. Matas, P. Priputen, K. Flachbart, K. 
Siemensmeyer, E. Wulf, A. Evdokimova, and N. Shitse- 
valova. Acta Phys. Pol. A 113, 227 (2008). 
K. Siemensmeyer, E. Wulf, H.-J. Mikeska, K. Flachbart, 
S. Gabani, S. Matas, P. Priputen, A. Evdokimova, and N. 
Shitsevalova, Phys. Rev. L ett. 101, 177201 (2008). 
Z. Y. Meng and S. Wessel, larXiv:0808.3T04] 
M. Levin and C. P. Nave, Phys. Rev. Lett. 99, 120601 
(2007). 

Y.-Y. Shi, L.-M. Duan, and G. Vidal, Phys. Rev. A 74, 
022320 (2006). 

M. Hinczewski and A. Nihat Berker, Phys. Rev. E 77, 
011104 (2008). 

Z. C. Gu, M. Levin, and X. G. Wen, Phys. Rev. B 78, 
205116 (2008). 

Other constructions have also been suggested in the lit- 
erature. For example, see F. Verstraete, M. M. Wolf, D. 
Perez-Garcia, and J. I. Cirac, Phys. Rev. Lett. 96, 220601 



(2006); R. Oriis and G. Vidal, Phys. Rev. B 78, 155117 

(2008). 



